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Abstract 



In this paper we des cribe how to imp r ove th e performance of the symbolic- 
numeric method in (|Li and Zhi . 2009 . 2011 ) for computing the multiplicity 
structure and refining approximate isolated singular solutions in the breadth 
one case. By introducing a parameterized and de flated system with smooth - 
ing parameters, we generalize the algorithm in (|Rump and GraillatL 120091 ) 
to compute verified error bounds such that a slightly perturbed polynomial 
system is guaranteed to have a breadth-one multiple root within the com- 
puted bounds. 
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1. Introduction 



It is a challe nge proble m to solve the polynomial systems with singular 
solutions. Rail (iRalll . Il96d ) studied the convergence properties of Newton's 
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in (Rump and Graillat. 


2009; 
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2011) 



they computed verified error bounds for singular solutions of nonlinear sys- 
tems. 

In (|Li and Zhil . l201ll ). we present a symbolic-numeric method to refine an 
approximate isolated singular solution of a polynomial system when the Ja- 
cobian matrix of the system evaluated at the singular solution has corank one 
approximately. Our approach is based on the regularized Newton iteration 
and the computation of differential conditions satisfied at the approximate 
singular solution. The size of matrices involved in our algorithm is bounded 
by the number of variables. The algorithm will converge quadratically if the 
approximate singular solution is close to the isolated exact singular solution. 
A preliminary implementation performs well in most cases. However, it may 
suffer from computing and storing dense multiplicity structures caused by 
linear transformation or dense ex pressi ons of differential functionals even for 
sparse input polynomials. In (Li, 201 ll ). we show briefly how to evaluate the 
differential conditions more efficiently by avoiding the linear transformation 
and solving a sequence of least squares problems. The techniques for con- 
structing a parameterized deflation system and evaluations of differe ntial 
conditions are similar to those introduced in (|Levkin et all l200d . l2007l ) . 



Main contribution. In this paper, we still focus on the special case where 
the Jacobian matrix has corank one. We describe how to preserve the sparse 
structure of the input polynomial systems by avoiding the linear transforma- 
tion. We reduce the storage space for computing the multiplicity structure 
by saving and evaluating differential conditions instead of explicit construc- 
tion of differential functionals. Fur therm ore, we show that the parameterized 
deflated system introduced in for avoiding the construction of the 

dif ferential functionals repea tedly can be used to generalize the algorithm 
in ( Rump and Graillat . 20091 ) to compute verified error bounds, therefore, a 
slightly perturbed polynomial system is guaranteed to have a breadth-one 
multiple root within the computed bounds. We prove that it is always pos- 
sible to construct a regular augmented system to compute an inclusion of 
the singular root by choosing properly smooth parameters and renumber- 
ing the polynomials. We provide numerical experiments to demonstrate the 
effectiveness of our method. 



Structure of the paper. Section[2]is devoted to recall some notations and 
well-known facts. In Section [3l we describe a new algorithm for computing 
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the multiplicity structure of the singular solution when the Jacobian matrix 
has corank one. Some experiment results are given to show the efficiency of 
the new algorithm. In Section 01 we show how to construct a parameterized 
deflated system to refine and compute verified error bounds for the breath- 
one multiple roots. Some numerical examples are given to demonstrate the 
performance of our algorithm. 



2. Notation and Preliminaries 



Let R = K[x] denote a polynomial ring over the field IK of characteristic 
zero. Let / = (/i, . . . , f n ) be an ideal of R, x G IK™ an isolated root of /, 
nT-x = (%l — x\, . . . , x n — x n ) the maximal ideal at x. Suppose Q x is the iso- 
lated primary component whose associate prime is m x , then the multiplicity 
H of x is defined as the dimension of the quotient ring R/Q x . 

Let : R — > IK denote the differential functional defined by 



l 



0H 



9 



-(*), V 5 (x)Gi?, 



1) 



ax\ •••«„! dx^---dxt 

for a point x G lK n and an array a G N n . The normalized differentials have 
a useful property: when x = 0, we have dg (x. 13 ) = 1 if a = (3 or otherwise. 
We may occasionally write d a = ci" 1 d% 2 ■ ■ ■ d% n instead of d? for simplicity 
if x is clear from the context, where df l 



1 d a i 
a t \ dx" 1 ' 



Definition 2.1. The local dual space of I at x is the subspace of elements 
o/23 x = Span K {d^,a G N n } that vanish on all the elements of I 

V x := {A G Dx | A(/) = 0, V/ G /}, (2) 

where dim(2? x ) = [i. 

Computing a closed basi s of the local dual space is done essentially by 

matr i x-kernel computations (IMarinari et all Il995l : iMourrainl . Il996l : iDayton and Zengl . 
20051 : IWu and Zhil . 120081 : IZend . l2009h . which are based on the stability prop- 
erty of P x : 



VAgP x , ^,(A)eI)! 



)H, i = l,...,n, (3) 

where T>t denotes the subspace of £> x of the degree less than or equal to t, 
for t G N, and $ Xi : S x —> S x are the linear anti- differentiation operators 
defined by 



-l,...,a n ) 



if a-i > 0, 
otherwise. 
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Lemma 2.2. llStetten , \2004 . Theorem 8.36) Suppose {Ai, . . . , A s } is a closed 
basis o/D^T 1 , then an element A 6 ID* lies in T>$ if and only if it satisfies 
(0j and A(fi) = for i = 1, . . . , n. 

In fact, ([3]) is equivalent to finding \ >k 6 K. such that A 6 S>x satisfying 

* ai (A) = Aj^Ai + A;, 2 A 2 + • • • + A 4jS A„ for % = 1, . . . ,n. (4) 



If A,; fc are known, we can compute A by the following formula (jMourrainl . 



1996) 



A = S A iJ*-i( A i) + S A 2,i^ 2 (Ai) + • • • + X] V^JA;), (5) 
j=i j=i j=i 



where the differentiation operators ^> Xi : 2? 

|(ai,...,Oi+l,...,a„) 




if ot\ = • • • 
otherwise. 



Si are defined by 
0. 



Here and hereafter, let Jf(x) denote the Jacobian matr ix of the poly- 
nomial system F evalu ated at x. It has been noticed in (jStanlevl . Il973l : 



Davton and Zend . 120051 ) that when the corank of the Jacobian matrix Jf(x) 



is one, T>% has the important property: 



dim(X>i) - dim(£> 



1, for 1 < t < n - 1. 



Hence, it is also called t he breadth one c ase in (jDavton and Zend . l2005l l. 
For this special case, in (|Li and Zhil . l2009h . under the assumption that the 
first column of Jf(x) is zero, we employ both normalization and reduction 
techniques to compute a closed basis of T>± very efficiently by solving jj, — 1 
linear systems with the size bounded by n x (n — 1). 

Theorem 2.3. Xn and Zh\ . 200f\ . Theorem 3.1) Suppose x is an isolated 
breath-one singular root of a given polynomial system F = . . . , f n } with 
the multiplicity \i, the first column of the Jacobian matrix Jf( x ) is zero. Set 
Ai = 1 and A 2 = d\, then we can construct A k incrementally for k from 3 
by 

Afc = A fc + a k ^d 2 + a k ^ 3 d 3 H h a kj7l d n , (6) 

where A k is a differential functional which has no free parameters and can 
be obtained from previous computed {Ai, A 2 , . . . ,A k -i} by 

k-l k-l 

i=2 3=2 
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The parameters a^j, for j = 2, . . . , n, are determined by solving 



J F (±) 



a>k,2 




' A & (/i) " 


0>k,n 




. Mfn) _ 



(8) 



where Jf (x) consists of the last n — 1 columns of Jf (x) . This process will 
be stopped if there is no solution for We get the multiplicity \i = k — 1 
and {Ai, A2, . . . , A-} a closed basis of the local dual space T>±. 

In ( Li and Zhi 20091 ). when the first column of Jf(x) is not zero, we 
apply a linear transformation of variables to obtain a new system and a 
new root, which will satisfy the assumptions of Theorem 12.31 Finally, we 
can derive a closed basis of the local dual space of the original system at 
the original root by transforming back the computed basis. Unfortunately, 
these transformations always result in dense systems even the original ones 
are sparse. 



EXAMPLE 2.1. lOjikl \l98i) Consider a polynomial syst 



em 



x\ + x 2 - 3, x\ + \x\ 



The system F has (1, 2) as a 3-fold isolated zero. 
The Jacobian matrix of F at (1,2) is 



J F (1,2) 



2 1 
1 i 



which has a non-trivial null vector r 
formation of the variables 

1 

X\ = — 
2 

to obtain a new dense polynomial system 



(— |, 1) T . We apply a linear trans- 
yi + 2y 2 ,x 2 =y\ + y 2 , 



G 



1 



1 



1 



-£Ji ~ 2ym + 4y 2 + yi + 2/2 - 3, -y x + -7/12/2 + g?/2 - 2 Vl + 22/2 ~ 2 ' 



The returned closed basis of the local dual space of G at the new point 
(f,f) T by Theorem O is 

1 



Ai = 1,A 2 



di,A 3 



d 



20 



:d 2 



5 



which can be transformed back to a closed basis of F at (1, 2) 



Ai = 1,A 2 



1 11 

-di + d 2 , A 3 = -d\ - -did 2 + d\ 



10 



20 



-do. 



3. A Modified Algorithm for Computing a Closed basis of the 
Local Dual Space 

In this section, we show how to avoid the linear transformations in com- 
puting a closed basis {Ai, A2, . . . , A^} of the local dual space V%. 

Let r = (ri,r 2 , ■ ■ ■ ,r n ) T be a non-trivial null vector of Jp(x). Without 
loss of generality, we assume 



l r i| > \ r j\, for 1 < j < n. 

Otherwise, one can perform changes of variables to guarantee 
Then we normalize r by r\ and derive that 



a 2 — [1, 0,2,2, • • • 1 0,2, r. 



1, 



n 



(9) 

is satisfied. 
(10) 



is also a non-trivial null vector of Jp(x). We set Ax = 1 and 

A 2 = d x + a 2 ,2d 2 H h a2,rA, 1 02,2 1 <!,-■■, 1 02.nl < 1- (H) 



Lemma 3.1. Under the assumption of (U\) , the differential functional mono- 
mial d±~ l does not vanish in A&, for k = 3, . . . , fi. 

Proof. If fj, > 3, A3 will satisfy ^ and (dJ), and at least one of Ax i 2,A 2 ,2 ! - • 
\ n ,2 is not zero. Let t be the integer such that Xt,2 7^ 0, then the differential 
functional monomial d\d± does not vanish in A3 according to (J5]). It follows 
that dt does not vanish in $ a;i (A 3 ). Then from @ and (jlip . we derive that 
Ax, 2 7^ 0. Therefore, (if does not vanish in A3. 

The rest proof is done by induction. Assume the lemma is true for k and 
k < fj,, then similar to the analysis above, A^+x satisfies ([3]) and ([4]). Let t be 
the integer such that At k 7^ 0, then d\~ l dt does not vanish in Aj, + x- It follows 
that d^~ 2 dt does not vanish in $> Xl (Ak + i). Since dim(Xtf) - dim(X>|- 1 ) = 1, 
degree(Afc) = /c — 1, we derive that Ax,fc 7^ 0. Therefore, d\ does not vanish 
inA fc+1 . □ 
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Remark 3.2. According to Lemma \3.1\ the coefficient of d^" 1 in A k is not 
zero, then can be normalized to be 1. Moreover, we can assume that A k 
does not have terms {1, d\,d\, . . . ,d^~ 2 } inside. Otherwise, one can reduce 
them by {Aj, A2, • • • , A^_i}. These normalization and reduction can help us 
reduce the number of free parameters in ^ to n — 1. 

Lemma 3.3. Under the assumption of and after performing the nor- 
malization and reduction above, we have 

$ Xl (A fc ) = A*_ a , _ 
^(Afc) = a kti Ai + a fc _i ] jA 2 H h a 2ii A fc _i, for i = 2, . . . ,n, 

where aj^ is the coefficient of d{ in Aj, for 2 < j < k and k < [i. 

Proof. By Lemma 13.11 and Remark 13. 2| we know that A*, has a term 
d^ -1 and there are no terms of {1, d\,d\, . . . , d±~ 2 } in A k . Hence, according 
to (dJ), we derive that <3? Xl (A fc ) = A fc _i. Furthermore, since ^(Ajfc) G £>~~ 2 , 
we have 

^(Afc) = Aj,iAi + A ij2 A 2 H h Aj )fe _iA fc _i, for i = 2, . . . ,n. 

Using ((5J), we claim that Ajj is equal to the coefficient of {d\ 1 di) in A&. On 
the other hand, we know that $|,7 (Afe) = Ak-j+i, hence the coefficient of 
(d\ di) in Afc is equal to the coefficient of di in Ak-j+i which is equal to 
a k-j+i,i- Hence, Ajj = ak-j+i,i for 2 < j < k — 1 and we prove the second 
equality in ([12]) . □ 



According to Lemma l3~3j from a closed basis {Ai, A 2 , . . . , A k -i} of D^~ 2 
to compute a new element A^ in /T>^ 2 , the only n — 1 free parameters 
are a^. Now we modify Theorem 12.31 under the assumption (fTT]) to avoid 
the linear transformations. Note that we also adopt a new equivalent form 
(|14p to construct instead of ([7]), which introduces an efficient method for 
evaluating the differential functionals. 



Theorem 3.4. |xj \20lA ) Suppose x is an isolated breath-one singular root 



of a given polynomial system F = . . . , f n } with the multiplicity /i. Set 
Ai = 1 and A 2 = d\ + a 2j2 o? 2 + • • • + a 2jn d n , f/ien we can construct A k 
incrementally for k from 3 07/ 

Afc = A fc + afc, 2 fi 2 + a fci 3d 3 H h a k>n d n , (13) 
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where 




The parameters akj, for j = 2, ... ,n are determined by solving |2j), When 
k = fi+1, there is no solution for fSjj and the process will be stopped. The 
set {Ai, A2, . . . , A^} is a closed basis of the local dual space T>^. 

Proof. According to Lemma [3T31 and (15|13|14p . the constructed A& satisfy 
the stability property ([3]). Moreover, solving flSJ) will guarantee that Afc(/j) = 
0, for i = l,...,n. Therefore, by Lemma 12.2^ the set {Ai , A2, • • • , A M } is a 
closed basis of the local dual space D*. □ 



Remark 3.5. As showed in (Li and Zhu 200d ), for solving only the 



vector on the right side is updated, while the matrix of the size n x (n — 1) 
on the left side is fixed. So we apply the LU decomposition to Jp(ic), then 
solve two triangular systems instead of {2JJ. 

Now we consider Example 12.11 again. Since Jf(1, 2) has a non-trivial null 
vector r = (—5, 1) T , we perform a change of variables x\ o x%, then apply 
the method described in Theorem 13.41 for computing a closed basis of the 
local dual space of F at (1,2). We derive that 

Ai = 1,A 2 = --di + d 2 ,A 3 = ~d\ - ^did 2 + d 2 2 - ~d x . 



In dLi and Zhil . Eo09), in order to compute A k of V^/V^' 2 , we need 



to construct by ([7]) and evaluate Afc(/j), for i = 1, . . . ,n. Even if the 
input system F is sparse, the differential functional could still be very 
dense. Hence, the evaluation of the vector on the right side of (|8|) could be 
very expensive sometimes. 

EXAMPLE 3.1. Consider a polynomial system F = {f\, . . . , f s } 
fi = xf + x 2 - x i+1 , ifi<s, 

fs = x s , 

with a breath-one singular zero (0, ... ,0) of the multiplicity 2 s . 



S 



As shown in (|Li and Zhil . l2009h . for s = 6, about 17MB of memory is 



used to store the local dual bases and it takes about 3 hours to compute all 
of them. Moreover, for s = 7, we are not able to obtain all in 2 days, 
and for s = 9, the estimated store space is about 1GB. It is not a surprise 
that the computation is dominated by the evaluation of A^(F) in ([5}. 

In fact, we can view and as diff erentiation operators, denoted by 
Pk and Lk respectively ( Li and Zhi . 20091 ) . then we can take advantage of 



(fTB"j) and dH]) to construct the polynomial systems Pk(F) and Lk(F) by 

fc-2 

P k( F ) = jpi ' J L k - 3 (F) ■ a i+ i and L k (F) = P k (F) + J F ■ a k , (15) 



where Jlj(f) 1S the Jacobian matrix of Lj(F), a 2 = [1,^2,2) • • • > a 2,n] T an d 
a.j = [0, dj t 2, • • • , aj,n] T for j = 3, . . . , k — 1. Hence, we can compute the eval- 
uation of Pk{F) which is equal to A^(F) without constructing and storing 
the dense differential functionals A&. 

The routine MSB1 below takes an ideal I = (/x, /2, ■ ■ ■ , f n ) C R and an 
isolated root x E K n of / as input, where the corank of the Jacobian matrix 
Jf(x) is one, and returns the multiplicity u and a closed basis of the local 
dual space Pj. Besides, we take a 2 , a3, . . . , as output too, since one can 
construct all by (fT3j) and (fl4l) immediately after they are computed. 
Another reason for outputting a^'s is that these values are important for 
multiple root refinement and verification if x is given with limited precision, 
which will be discussed in the next section. 

Algorithm 3.6. MSB1 

Input: A polynomial system F = {fi, / 2 , . . . , /„} and a root x E K n . 
Output: The multiplicity p, the parameters a2, a3, . . . , a„ and a closed basis 
{Ai, A2, . . . , A^} of the local dual space T>%. 

1. Compute a non-trivial null vector r = [r%, r 2 , ■ ■ ■ , r n ] T of Jp(x). Let t 
be the integer s.t. \rt\ > \rj\, 1 < j < n. 
Apply variables exchange x% -H- x% to F, x, Jp and r. Set 



a 2 := 



, r2 rs rn 

? ? ' ' ' 3 

ri r\ r\ 



, L 2 (F) := Jp ■ a 2 and P 3 := ^Jl 2 (F) • a 2 . 



Compute the LU Decomposition of Jp(x) = P ■ L ■ U . Set k := 3. 
2. Compute Pk(F) by [W\) and evaluate it at x io get Afc(F), and sofee 

L • b = -P- 1 • A fc (F). 
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If the last entry in b is zero, solve f7i..( n -i), : ■ c = b, and set 

and L k (F) := P k {F) + J F ■ a fc , 




c 



and repeat with k := k + 1. Otherwise, set [i := k — 1, go to 5top0 



3. Construct {Ai, A 2 , . . . , A^} fry ((J^j and ( fm ) using a 2 ,a 3 , 
crementally, then perform variables exchange d\ ^ d t . 



■ ■ ■ i a /j 



It should be noticed that we run Step [3] in Algorithm 13.61 only when 
a closed basis for the local dual space is wanted. In general, we can omit 
this step in refining and certifying approximate singular solutions. We have 
implemented this method in Maple. In the following table, we show the time 
needed for computing all a,- in Example 13.11 for j = 2, . . . , 2 s . 



s 


6 


7 


8 


9 


10 


multiplicity 


64 


128 


256 


512 


1024 


time(sec) 


0.593 


1.377 


3.445 


10.913 


44.659 



4. Verified Multiple Roots of Polynomial Systems 



As mentioned in (jMantzaflaris and Mourrainl . 120111 ) , in real-life applica- 
tions it is common to work with approximate inputs, and usually we need to 
decide numerically whether an approximate system possesses a unique real 
or complex root in a given domain. 

Standard verifi cation methods for nonl i near systems are b ased on the 
following theorem ( Krawczyk . 1969 : Moore . 1977 : Rump . 19831 ). 



Theorem 4.1. Let F G R n be a polynomial system with F = (fx, . . . , f n ), 
and x E ]R n a real point. Given an interval domain X € IK n with x £ X, 
and an interval matrix M £ IR nxn satisfies V/j(X) C Mi : , for i = 1, . . . , n. 
Denote by I the n x n identity matrix and assume 

- Jp(x)F(x) + (J - J F (x)M)X C int(X). 

Then there is a unique x£l with -F(x) = 0. Moreover, every matrix M £ 
M is nonsingular. In particular, the Jacobian matrix Jp(x) is nonsingular. 



In ([Rump and GraillatJ . l2009h . they introduced a smoothing parameter 
to certify a double root of a slightly perturbed system using Theorem 14.11 
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It should be noticed that a double root is the simplest breath-one root with 
the mu ltiplicity 2. 



In (jMantzaflaris and Mourrainl . l201ll ). they applied Theorem 14.11 to a 



deflated system to verify a multiple root of a nearby system with the com- 
puted local dual structure. Their method can deal with arbitrary multiple 
roots. 

For the breath one case and /i > 2, in ( Rump and Graillat . 20091 . Theo- 
rem 4.2), they proved that it is impossible to compute an inclusion of a mul- 
tiple root by adding only a smoothing parameter to one selected equation. 
We show below how to construct a deflated system using the parameter- 
ized basis in a2, . . . , a„ for the local dual space of T>± to certify breath-one 
multiple roots for fi > 2. 

Let F = {fx, . . . , f n } 6 R be given. Suppose x 6 K. n is an isolated root 
of F with the multiplicity fi and Jp(x) has corank one. We show first how to 
choose a pair of suitable variable and equation to perform the perturbation. 
In fact, as we showed before, the perturbed variable Xi can be determined by 
choosing a column of Jp(x), which can be written as a linear combination 
of the other n — 1 columns. Similarly, suppose the j-th row of Jp(x) can 
be written as a linear combination of the other n — 1 linearly independent 
rows, then we add the perturbed univariate polynomial in Xi to fj. Then 
we perform 

xx -H- Xi and fx O fj (16) 
to construct the deflated system in (fTT)) , 

Assumption 4.2. Suppose Jp(x) has corank one. We assume below that 
the first row (column) of Jp(x) can be written as a linear combination of its 
other rows (columns). This can always be achieved by changing of variables 
and renumbering equations as above. 

We introduce fJ>—l smoothing parameters bo,bx, ■ ■ ■ , and construct 
a deflated system G(x, b, a) with fin variables and fin equations: 



G(x,b,a) 



iq(x,b) 



= ^(x)-(E^ 25 #Jei 

F 2 (x,b,a 2 ) 
F 3 (x,b,a 2 ,a 3 ) 



V 



F M (x,b,a 2 , . . . ,a p 



(17) 



/ 



where b = [b , b x , ■ ■ ■ , 6^-2], a = [a 2 , a 3 , . . . , a M ], a 2 = [1, a 2)2 , . . . , a 2 , n } T , 
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afc = [0, a k ,2, ■ ■ ■ , a k , n ] T for 2 < k < fi, and 

F fe (x,b,a 2! ...,a fe ) =L k {Fx). (18) 



Theorem 4.3. Suppose G(x,b,a) = 0. Under Assumption if the Ja- 
cobian matrix J<j(x, b,a) is nonsingular, then x is an isolated root of the 
polynomial Fq(x) = Fi(x, b) with the multiplicity ^ and the corank of Jp (x) 
is one. 

Proof. From G(x, b,a) = 0, we have Fq(-x) = and 

F 2 (x, b, a 2 ) = Jf (x) • a 2 = 0. 

Since a 2 ^ 0, we derive that 

rank(J^ (x)) < n — 1. 

Moreover, from the expression of a 2 , we know that the first column of Jf (x) 
can be written as a linear combination of the other n — 1 columns. Therefore, 

rank(J^ (x)) = rank(J^ (x)) < n — 1, 

where Jf (x) consists of the last n — 1 columns of Jf (x) . Similarly, since 
Ffc(x, b, a 2 , . . . , afc) = 0, by Theorem 13.41 and (fT3)) . we derive that 

rank(Afc(Fo), J^ (x)) = rank(J^ (x)) < n — 1, for 2 < k < fi. 

In order to prove that x is a breadth-one root of Fq(x) = with the multi- 
plicity fi, we need to show that 

rank( Ji? (x)) = n — 1 and rank(A At+ i(i ? o), Jf (x)) = n. (19) 

It is interesting to see that we can use the equivalent relations 

()hl - ° l for 1 < i < k and 1 < j < n, (20) 



daij dxj 



to obtain a simplified expression of Jq 

( It, ... 2i p-, n n ... n n \ 



G 



J Fl e x xiex 
Jf 2 ei 



-e-i 














-ei 













-ei 


Jf 2 











"2 u e i JJTW) 
J F 3 ••• ^ 



J> M _i ••• ei Jf m _ 2 «/f m _ 3 ••• «[f 
\ Jf^ • • • Jp M _ 1 Jf m _ 2 • • • Jf 2 / 



(21) 
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where Jp k denotes the Jacobian matrix of Ffe(x, b,a2, . . . , a^) with respect 
to x and Jp. consists of the last n — 1 columns of Jp k , for 2 < k < [i. 

If rank( Jp (x)) < n — 2, there will exist a nontrivial vector in its kernel. 
Note that J^ (x) is the only non-zero element in the last column, we claim 
that there will exist a nontrivial vector in the kernel of Jg(x, b, a), which is 
a contradiction. Then we derive that 

rank(J Fo (x)) = rank(A fc (F ), Jf (x)) = rank(Jp (x)) = n - 1. 

On other hand, from f|14j> and (|18|) . we derive that 

J G (x,b,a). jL (/ ,_ 1)n+1 • v = (0, . . . ,0, A^ + i(F )) T , 

where 

v = — (1, a 2 ,2, «2,n, 0, . . . , 0, 2a 3j2 , . . . , 2a 3i „, . . . , 1)6 M ,2, • • • , (/x-l)a Atjn ) T . 
A 4 

So that, if rank(A At+1 (i ? ), Jf (x)) < n — 1, there will exist a nontrivial 
vector in the kernel of Jg(x, b, a), which is a contradiction. Hence, we have 

rank(A^ + i(F ), Jf (x)) = n. 

Therefore, according to Theorem 13.41 x is an isolated breath-one singular 
root of -Fo(x) = with the multiplicity (i. 

□ 

Theorem 4.4. Suppose x is an exact isolated root of F(x) = with the 
multiplicity fi and the corank of J_f(x) is one exactly. Under Assumption 
\4-S\ we have 

rank( Jp(x), ei) = n. (22) 

Evaluating ji7| ) at b = 0, i.e., no perturbations for F , the Jacobian matrix 
Jg(x, 0,a) vis nonsingular. 

Proof. According to Assumption 14.21 the first row of Jf( x ) can be written 
as a linear combination of its other rows. Since the rank of Jf(x) is n — 1, 
its last n — 1 rows must be linear independent. Therefore, we have (|22p . 

Assume v is a nontrivial vector in the kernel of Jg(x, 0,a). If v\ = 0, 
by checking the columns of Jq in ([2ip . and using (|22p . we can show that 
v = 0. If v\ 7^ 0, we can assume v\ = 1. Similar to the second part of 
proof of Theorem 14.31 we derive that A^ + i(F) can be written as a linear 



13 



combination of the columns from Jp(x), which is a contradiction. Hence, 
there exists no nontrivial vector in the kernel of Jg(x, 0, a). In other word, 
Jg(x, 0,a) is nonsingular. □ 

Now, we apply Theorem 14. II on the deflated system. If the test succeeds, 
we will derive verified and narrowed error bounds with the property that 
a slightly perturbed system is proved to have a breadth-one multiple root 
within the computed bounds. 



Theorem 4.5. Suppose Theorem \4.1\ is applicable to G(x, b, a) in {11) and 
yields inclusions for x, b and a such that G(x, b, a) = 0. Then x is an 
isolated breath-one root o/Fq(x) := Fi(x, b) with the multiplicity [x. 



Proof. A direct result of Theorem 14.11 and Theorem 14.31 



□ 



EXAMPLE 4.1. \Rump and Graillai 1 20091 . Example 4. 11) Consider a 
polynomial system 

7-1 r 2 2 2 1 

t = \X l X2 — XlX 2 , X\ — X 2 - j 

The system F has (0, 0) as a 4-fold isolated zero. 



The Jacobian matrix of F at (0, 0) is 

MQ,o) = 





1 



so we choose x 2 as the perturbed variable and add the univariate polynomial 
60 + bix 2 + y x 2 t° t ne first equation in F to construct the parameterized 
deflated system 



G(x,b,a) 



x\x 2 — x\x\ — h 



'0 



bix 2 - ^x\ 



\ 



\ 



2a\X\x 2 — a\x\ + x\ — 2x\x 2 — b\ — b 2 x 2 
a\ — 2x 2 

a\x 2 + 2a\X\ — 2a\x 2 + 2a 2 x\x 2 — a 2 x\ — x\ — -# 

0,2 — 1 

a\ + 2a\a 2 x 2 — a\ + 2a 2 x\ — 2a 2 x 2 + 203X1X2 — 03X2 

03 



Applying the INTLAB function verify nlss ( Rump . 19991 ) to G with the initial 
approximation 



[0.002, 0.003, -0.001, 0.0015, -0.002, 0.002, 1.001, -0.01] 
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to obtain inclusions 



-0.00000000000001 , 0.00000000000001] 

-0.00000000000001, 0.00000000000001] 

-0.00000000000001 , 0.00000000000001] 

-0.00000000000001, 0.00000000000001] 

-0.00000000000001, 0.00000000000001] 

This proves that the perturbed system -Fo(x) (\bi\ < 10 -14 , i = 0, 1, 2) has a 
4-fold root x with -10" 14 < x < 10" 14 . 

EXAMPLE 4.2. \Rump and Graillai\200<\ . Example 4-7) Consider a poly- 
nomial system 

17 r 2 2 2 l 

F = \X\ — x 2 ,xi - x 2 .\ 
The system F has (0, 0) as a 2-fold isolated zero. 

For this example, as mentioned in (|Rump and Graillatl . l2009h . the iteration 
is sensitive to the initial approximations. Applying the INTLAB function 
verifynlss2 to F with the starting point [0.002,0.001], we will obtain inclu- 
sions 



[ -0.00000000000001, 0.00000000000001] 

[ -0.00000000000001, 0.00000000000001] 

[ -0.00000000000001, 0.00000000000001] 

However, for the initial approximation [0.001,0.001], we obtain 



0.49999999999999, 
0.70710678118654, 
-0.25000000000001, 



0.50000000000001] 
0.70710678118655] 
-0.24999999999999] 



which finds the double root (0.5, 1/V2) of x\ — x| + 0.25 = and x\ — x\ = 0. 

For this reason, w e prefer to use the symbolic-numeric method described 
m (|Li and Zhil . l201lh to refine initial approximations firstly, then use the 
method in Theorem l4.5l to compute inclusions of multiple roots. We show the 
routine MRRB1 below for refining a singular solution to high precision. The 
input of MRRB1 is a sequence of polynomial systems F\, F 2 , . . . , F^ defined 



m 



W> 03 with b = 0, F, 



.P/,,-1-1 (F ) an d an approximate solution of 



Fi = 0. The algorithm in (|Li and Zhil . l201ll l has been improved in MRRB1 
by avoiding linear transformations and constructing differential functionals 
repeatedly. 
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Algorithm 4.6. MRRB1 

Input: A sequence of systems F\, . . . , Fu+i, a point x £ ~K n . 
Output: A refined point x and refined parameters a 2 , ■ ■ ■ , a„ . 

1. Regularized Newton Iteration: Solve the least squares problem 

(J^(x) • J Fl (x) + a n I n ) y = -J* Fx {±) ■ Fi(x), 

where J Fi (x) is the conjugate transpose of Jf 1 (x), a n is the smallest 
singular value of Jf 1 (x) and I n is the n x n identity matrix. 
Set x := x + y. 

2. For 2 < k < jj,, solve the least squares problem 



F fc (x,a 2 , . . . ,a fc _i,a fe ) 







to obtain a^. 
3. Solve the linear system 



F M+1 (x,a 2 , . 



<9Fi(x) 



5 <W> 



dF^yi) 



-^(x,a 2 , 



where v = (v\, . . . , v n ) T . Set 5 := v 
4. Return a 2 , . . . , a„ and 



x := x + 5 



( 1 \ 

02,2 



\ "2,n / 



Now we consider Example 14.21 again. For [0.002, 0.001], after running MRRB1 
two times then applying the INTLAB function verify nlss to G, we obtain 

[ -0.00000000000001, 0.00000000000001] 
[ -0.00000000000001, 0.00000000000001] 
[ -0.00000000000001, 0.00000000000001] 

Similarly, for [0.001,0.001], we obtain 

[ -0.00000000000001, 0.00000000000001] 
[ -0.00000000000001, 0.00000000000001] 
[ -0.00000000000001, 0.00000000000001] 
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Xi + Xi - x i+ i, if i < s, 

„3 



EXAMPLE 4.3. Xhi and Zhi \200d) Consider a system F = {/i, . . . , /,} 
given by , 

Ji 

fs = x„ 

with a breath-one singular zero (0, 0, ... , 0) of the multiplicity 3. 

We run MRRB1 three times for initial approximate roots near the origin, 
whose errors are around 10 -4 , to obtain refined x and a with errors about 
1CP 12 . We choose x s as the perturbed variable and add the univariate poly- 
nomial &o + b±x s to the last polynomial f s to construct the parameterized 
deflated system. In the following table, \X\ and \B\ denotes the interval size 
of inclusions for x and b computed by applying INTLAB function verify nlss 
to the deflated system (fT7|) and (x, 0,a). 



s 


\x\ 


\B\ 


10 


10 -14 


10 -14 


20 


10 -14 


10 -14 


50 


io- 14 


10 -14 


100 


10 -14 


1Q -14 


200 


10" iy 


io- 112 


500 


10" iy 


io- iy 


1000 


io- 12 


io- 12 
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